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[i]  Techniques  are  presented  for  efficiently  evaluating  quasi-linear  diffusion  coefficients 
for  whistler  mode  waves  propagating  according  to  the  full  cold  plasma  index  of  refraction. 

In  particular,  the  density  ratio  upe/Q.e  can  be  small,  which  favors  energy  diffusion.  This 
generalizes  an  approach,  previously  used  for  high-density  hiss  and  electromagnetic  ion 
cyclotron  waves,  of  identifying  (and  omitting)  ranges  of  wavenormal  angle  8  that  are 
incompatible  with  cyclotron  resonant  frequencies  to  occurring  between  sharp  cutoffs  of  the 
modeled  wave  frequency  spectrum.  This  requires  a  detailed  analysis  of  the  maximum 
and  minimum  values  of  the  refractive  index  as  a  function  of  to  and  8,  as  has  previously 
been  performed  in  the  high-density  approximation.  Sample  calculations  show  the  effect  of 
low-density  ratio  on  the  pitch  angle  and  energy  diffusion  coefficients  modeling  the 
effect  of  chorus  waves  on  radiation  belt  electrons.  The  high-density  approximation  turns 
out  to  be  quite  robust,  especially  when  the  upper  frequency  cutoff  is  small  compared  with 
ffy  The  techniques  greatly  reduce  the  amount  of  computation  needed  for  a  sample 
calculation,  while  taking  into  account  all  resonant  harmonic  numbers  n  up  to  ±oo. 

Citation;  Albert,  J.  M.  (2005),  Evaluation  of  quasi-linear  diffusion  coefficients  for  whistler  mode  waves  in  a  plasma  with  arbitrary 
density  ratio,  J.  Geophys.  Res.,  110,  A03218,  doi:10.1029/2004JA010844. 


1.  Introduction 

[2]  Cyclotron-resonant  wave-particle  interactions  are  cur¬ 
rently  a  leading  candidate  for  understanding  the  behavior  of 
outer  radiation  belt  electrons  following  magnetic  storms.  In 
particular,  Summers  et  al.  [1998]  proposed  that  energization 
by  dawnside  whistler  mode  chorus  waves,  combined  with 
pitch  angle  scattering  by  duskside  electromagnetic  ion 
cyclotron  (EMIC)  waves,  might  account  for  the  observed 
distributions  of  energetic  electrons.  By  analyzing  the  path  in 
velocity  space  followed  by  an  electron  as  it  maintains 
primary  (n  =  —1)  resonance  with  a  field-aligned  wave,  it 
was  found  that  energy  diffusion  is  favored  by  the  low  (~1) 
values  of  the  density  ratio,  wpe/$le,  found  outside 
the  plasmasphere.  Home  et  al.  [2003]  presented  detailed 
quasi-linear  diffusion  coefficients  for  low-density  chorus 
and  confirmed  the  effectiveness  of  energy  diffusion,  al¬ 
though  these  calculations  were  restricted  to  local  (equatorial) 
interactions.  Horne  et  al.  [2005]  presented  bounce-averaged 
results,  using  realistic  wave  models.  These  calculations 
were  performed  using  two  different  approaches,  one  of 
which  is  described  by  Glauert  and  Home  [2005]  and  one 
of  which  is  described  here. 

[3]  Whistler  mode  chorus  waves  are  observed  to  have 
fine  structure,  including  time-evolving  frequencies  [e.g., 
Santolik  et  al.,  2003]  and  are  associated  with  the  strongly 
nonlinear  processes  of  particle  phase  bunching  and  phase 
trapping  (see  the  review  by  Sazhin  and  Hayakawa  [1992], 
and  also  Smith  and  Nunn  [1998],  Albert  [2000,  2001, 2002], 
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and  Trakhtengerts  et  al.  [2003]).  However,  it  is  assumed 
here  that  quasi-linear  theory  is  a  meaningful  description  of 
the  long-term  effect  of  fully  formed  waves  on  high-energy 
particles  [Summers  et  al.,  1998;  Summers  and  Ma,  2000; 
Home  et  al,  2003,  2005]. 

[4]  Previous  calculations  of  quasi-linear  diffusion  have 
generally  relied  on  the  high-density  approximation  intro¬ 
duced  by  Lyons  [1974b],  which  leads  to  a  greatly  simplified 
form  of  the  index  of  refraction  p.  This  allows  a  character¬ 
ization  of  the  resonant  frequencies  at  each  wavenormal 
angle.  Since  the  diffusion  coefficients  are  infinite  sums  of 
integrals  over  wavenormal  angle  8  (and  then  bounce  aver¬ 
aged),  it  turns  out  to  be  possible,  and  extremely  beneficial, 
to  restrict  the  range  of  the  integrals  over  8.  This  analysis  has 
been  carried  out  for  whistler  mode  plasmaspheric  hiss 
[Albert,  1999,  hereinafter  referred  to  as  Paper  1]  and  for 
EMIC  waves  [Albert,  2003,  hereinafter  referred  to  as  Paper 
2],  using  the  high-density  approximation.  This  paper  gen¬ 
eralizes  that  approach  to  whistler  mode  waves  with  arbitrary 
values  of  the  density  ratio. 

[5]  Section  2  presents  the  quasi-linear  diffusion  coeffi¬ 
cients  (pitch  angle,  momentum,  and  mixed),  in  a  form  that 
generalizes  the  formulas  of  Lyons  [1974b]  to  relativistic  test 
particles,  general  p,  and  arbitrary  wavenormal  angle  and 
frequency  profiles.  Section  3  expresses  the  cyclotron  reso¬ 
nance  condition  in  terms  of  1/p2  (denoted  by  4/)  and  a 
function  V  which  contains  all  the  dependence  on  parameters 
of  the  test  particle  and  briefly  discusses  the  simple  geomet¬ 
ric  behavior  of  V  as  a  function  of  w.  Section  4  presents  a 
detailed  algebraic,  geometrical,  and  numerical  analysis  of 
the  full  two-component  cold  plasma  whistler-mode  refrac¬ 
tive  index  as  a  function  of  w.  Particular  attention  is  paid  to 
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the  maximum  and  minimum  values  over  a  fixed  w  range  and 
how  these  values  vary  with  9.  Several  bounds  and  approx¬ 
imations  are  discussed,  including  the  high-density  approx¬ 
imation  and  a  novel  factorization.  Section  5  identifies  values 
of  9  for  which  the  maximum-to-minimum  ranges  of  F(w) 
and  'f'(oj)  do  not  overlap;  for  such  9  values,  there  can  be  no 
resonances.  For  large-enough  harmonic  number  |n|  there 
will  be  no  resonances  for  any  9  (since  Vmax  becomes  less 
than  ^min),  which  cuts  off  the  sum  over  n  =  — oo  to  oo. 
Additional  tests  based  on  the  curvatures  of  F(w)  and  'I'(w) 
are  also  developed.  Finally,  the  dependence  on  latitude  X  is 
considered,  which  leads  to  tests  for  X  ranges  that  have  no 
resonances  for  any  9  or  n.  In  section  6  the  preceding 
analysis  is  used  for  n  <  0  to  isolate,  or  bracket,  the  resonant 
frequencies  which  do  exist  so  that  they  can  be  found  quickly 
and  easily  by  a  generic  one-dimensional  (1-D)  root- finding 
algorithm.  This  bracketing  cannot  always  be  done  for  n  >  0, 
as  explained,  and  alternative  procedures  are  discussed. 
Section  7  discusses  the  ratio  Dpp/Daa  and  explores  why  this 
ratio  should  increase  as  upJQ.e  decreases.  Section  8  explores 
the  reliability  of  the  high-density  approximation  for  different 
values  of  upe/ne  and  quantifies  the  gain  in  computational 
efficiency  made  possible  by  the  methods  of  this  paper.  A 
summary  is  given  in  section  9. 

2.  Quasi-Linear  Diffusion  Coefficients 

[6]  The  resonance  condition  for  a  diffusing  particle  is 

U)  -  *|| V||  =U„  (I) 

where  n  is  an  integer,  s  =  ±1  is  the  sign  of  the  charge  of  the 
particle,  fic  =  \q\BJmc  is  its  local  nonrelativistic  gyrofre- 
quency,  and  y  is  its  relativistic  factor.  The  local  pitch  angle 
of  the  particle  is  a,  the  index  of  refraction  is  pi  =  kcl w,  the 
wavenormal  angle  is  9,  and  x  will  denote  tan  9  throughout. 

[7]  The  underlying  physics  of  the  diffusion  coefficients 
used  here  is  identical  to  that  of  Lyons  [1974b],  except  for 
the  extension  to  relativistic  test  particles  (as  briefly  dis¬ 
cussed  by  Lyons  [1974a]).  However,  the  frequency  and 
wavenormal  angle  profiles  ZI2(w)  and  gw(9)  are  kept  general, 
and  we  also  postpone  specifying  the  form  of  the  refractive 
index.  The  normalized  diffusion  kernel  for  Daa  can  be 
written  as 


The  normalized  diffusion  coefficients  themselves  are  then 


®-=e/^>=e/^ 


dx 


£/u,  (5) 


where  n  ranges  over  ±00,  and  similarly  for  Vrjp  and  Vpp. 
The  term  above  is  exactly  as  given  in  equation  (9)  of 
Lyons  [1974b],  except  that  the  arguments  of  the  Bessel 
functions  therein  become  [nxtan  ct(w  —  w„)/w„]  (as  in  Paper  1 ) 
to  account  for  relativistic  factors.  The  partial  derivatives 
above  may  be  expressed  in  terms  of  derivatives  of  p(w,  9)  by 


and 


1  / _  w  p. 

v||  w  -  w„  p  4-  ud\i/du 


(6) 


dx  _  pw„/[u(w  —  u>„)]  -  d\i/du) 
3w  psin8cos0  -  d\i/dx 


(7) 


The  physical  momentum  diffusion  coefficients  D  are  equal  to 
the  normalized  quantities  V  multiplied  by  nc(5^ave/5o)(/nv)2. 
Finally,  the  bounce-averaged  diffusion  coefficients  are 
=  /  Daa(dao/d(x)2  dthb,  (D^p)  =  f  D^d a^da) 
dthb ,  and  (Dpp)  =  J  Dppdihb,  where  rb  is  the  bounce  time. 
These  may  be  converted  to  integrals  over  latitude  (as  detailed 
for  Daa  by  Lyons  et  al.  [1972]  for  a  dipole  field)  or  over 
a  (which  may  give  better  numerical  behavior  near  the 
mirror  point). 

fa]  It  has  become  standard  to  model  the  power  spectrum 
2r( w)  as  proportional  to  a  gaussian  truncated  at  u >£,c  and 
w uc  and  the  wavenormal  angle  profile  gu(x)  as  a  gaussian 
truncated  at  and  xmax.  Then,  taking  Aw  equal  to  the 
frequency  width  6w, 


G,  = 


2  exp[-(w-  /bJ 

vAt  erf  [(wye  -  wn)/Bw]  +  erf[(wm  -  wiC)/6w]  ’ 


(8) 


and,  with  the  high-density  approximation  for  p  (discussed 
below), 


/  M  3/2 exP  -(*  ~  xmf/xj 

i{  w) 


_  irnc  cos2e  ^  (-  sin2  a  +  w„/u)2 
aa  2Aiv„/c|3  n|i-(^/a*||)I/v|||  1  21 


where  G,(w)  =  (Aw)S2(w)/ / B2{J)dJ,  and 

gu,(8) 


G2(w,0)  = 


J  gu(0,)p2|p+w9p/9w|  sinOW 


(3) 


The  normalization  Aw  is  arbitrary,  but  it  is  convenient  to  use 
a  value  characteristic  of  the  width  of  32(w).  Then  V^a, 
and  Vpp  are  given  by 


sin  a  cos  a  (V^V 

—  sin2  a  +  w„/u)’  UJ  ' 


where  M  =  mjmi  and  /(w)  (with  xmin  =  0)  is  defined  by 
Lyons  [1974b].  Nonrelativistically,  these  expressions  all 
reduce  exactly  to  the  formulas  given  there.  The  choice  of 
gaussians  is  not  crucial,  but  the  existence  of  lower  and 
upper  frequency  cutoffs,  wic  and  wuc,  is  essential  for  the 
analysis.  As  in  Papers  1  and  2,  the  frequency  cutoffs  can  be 
used  to  derive  restrictions  on  9,  which  limit  the  range  of  the 
x  integration  in  equation  (5).  Thus  by  preidentifying  ranges 
of  9  for  which  w  is  outside  the  cutoffs,  the  amount  of 
computation  needed  to  obtain  D  may  be  drastically  reduced. 

3.  Resonance  Condition 

[9]  Despite  the  formidable  appearance  of  the  expressions 
for  the  diffusion  coefficients,  the  only  quantity  not  straight¬ 
forward  to  evaluate  is  the  resonant  frequency  w(9),  which  is 
needed  if  V**  is  integrated  over  9  rather  than  w. 
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V(w) 


u/ne 

Figure  1.  The  function  F(uj)  at  fixed  0,  with  n  =  0  (dotted 
curve),  n  =  1  (dashed  curve),  and  n  =  —  1  (solid  curve),  for 
an  electron  with  7  =  4  and  (vf/c2)cos2  0  =  0.5. 

[10]  The  resonance  condition,  equation  (1),  may  be  writ¬ 
ten  as  V(bj,  0)  =  ^(ui,  0)  where,  similar  to  Papers  1  and  2, 

F=-^ — — — r  cos2 acos2 0,  'J'  =~.  (10) 

c2  (w  -  w„)2  P2 

There  is  a  loss  of  sign  information  in  going  from  equation 
(1)  to  V=  \P  (via  squaring),  but  it  is  assumed  that  for  every 
wave  mode  with  positive  k\\/k  there  is  a  mode  of  equal 
amplitude  with  negative  k^/k,  so  any  solution  of  V  =  \P 
corresponds  to  a  solution  of  (1).  Furthermore,  the  particle 
distribution  is  assumed  to  be  the  same  for  vjj  positive  and  V|| 
negative.  This  definition  of  'I'  differs  slightly  from  that  of 
Lyons  [1974b]  and  Papers  1  and  2,  where  it  denoted  the 
high-density  limit  of  [(Wp<7fi2)(l  +  M)/M]/p.2.  This  limit  will 
be  discussed  in  the  next  section. 

[11]  The  geometric  properties  of  F(oj)  at  fixed  0  were 
established  in  Paper  1  and  are  shown  in  Figure  1 .  For  n  =  0, 
F(u>)  has  the  constant  value  F0  =  (v|/cr)cos2  0.  For 
positive,  F  always  has  positive  curvature  and  is  increasing 
for  uj  <  and  decreasing  and  bounded  below  by  F0  for  uj  > 
w„.  For  negative,  F(w)  is  strictly  increasing  but  bounded 
above  by  F0  and  has  positive  curvature  d2  V/dw2  for  w  < 
|uj„|/2  and  negative  curvature  for  w  >  |uj„|/2.  The  properties 
of  $  are  discussed  in  the  next  section. 


4.  Index  of  Refraction 

[12]  The  standard  wave  coefficients  [Stix,  1962]  for  a  cold 
pure  electron-proton  plasma  are 


P=  1 


=  1  ± 


uj2 


1  +M 


it, 

.,2 


fieU)  1  —  M  (u/Ue  —  fl./uj)  ’ 
R-L 


D-- 


and  the  full  expression  for  'I'  can  be  written  as 
$  =  I  ( RL-PS )  sin2  8  +  2 PS 


:  ^ / (RL-PS)2  sin4  0  +  4 P2D2  cos2  0 1/2 PRL. 


(11) 


(12) 


The  choice  ±  in  equation  (12)  determines  the  labeling  of  the 
wave  as  “R  mode”  (+)  or  “L  mode”  (— ),  but  the 
wavenormal  surfaces  exchange  labels  (or,  as  Stix  [1962] 
puts  it,  the  labels  exchange  wavenormal  surfaces)  when 
crossing  P  =  0  or  D  =  0.  Here  only  waves  with  w  <  and  uj 
<  u>pe  are  considered,  although  upe  is  allowed  to  be  less  than 
n.  Thus  P  is  always  negative,  but  D  becomes  negative  for 
uj  <  fi/  (to  lowest  order  in  M).  Therefore  ±  is  chosen  to  be 
“+”  when  D  >  0,  for  R-mode  polarization,  and  ”  when 
D  <  0,  for  continuity. 

[13]  Usually,  p,2  is  treated  as  a  function  of  0  at  fixed  uj; 
here,  Figure  2  shows  versus  uj  at  fixed  0  (and  fixed  ijjpJQ.e). 
For  each  curve,  it  can  be  shown  algebraically  that  d^/cluj  =  0 
at  uj  =  0  so  that  the  behavior  of  $  near  zero  is 
determined  by  its  curvature  d^'i/d uj2  there.  For  small 
and  moderate  0  this  curvature  is  positive,  so  'P(io)  starts 
out  increasing  but  soon  develops  negative  curvature, 
eventually  peaks  (before  reaching  He/2),  and  then 
decreases  to  zero  while  retaining  negative  curvature!  These 
curves  each  have  one  internal  maximum  and  no  internal 
minima.  For  large  0,  'P  starts  out  with  with  negative 
curvature  and  declines  monotonically  •  to  zero,  with  no 
internal  maxima.  The  two  classes  of  curves  are  separated 
by  the  curve  that  has  zero  curvature  at  w  =  0;  this  occurs  for 


sin2  0, 


'sep  ■ 


2(Rj 

RR" 


u— 0 


M  +  uj IJQ] 


|2  » 


(13) 


which  gives  0  very  close  to  90°  unless  upe/Qe  ~  \/M  or  less. 
This  is  shown  as  the  thick  curve  in  Figures  2c  and  2d  and  is 
observed  to  cross  zero  at  frequency  ujscp  near  Q.gm  = 
(SlSlj)'12  (the  geometric  mean  of  and  O,). 

[14]  In  general,  $  reaches  zero  at  the  frequency  given  by 
the  expression  tan2  0  =  —P/S  [Stix,  1 962].  Using  this  to  find 
the  $  curve  that  reaches  zero  right  at  for  which  5=1, 
yields 


sin2 


-P 
1  -P 


Qgm 


(14) 


which  is  slightly  less  than  sin2  0sep.  This  is  shown  by  the 
upper  dashed  curve  in  Figure  2d.  The  $  curve  for  sin2  0  =  1 
crosses  zero  at  the  lower  hybrid  resonance,  wLH,  where  S  = 
0.  Thus  for  w  >  u )Lh  there  is  some  maximum  0,  the  well- 
known  resonance  cone  angle,  beyond  which  $  becomes 
negative.  This  $  curve  is  the  lower  dashed  curve  in  Figure  2d. 
Since  sin2  0gm  <  sin2  0scp  <  1,  the  separating  curve  crosses 
zero  at  with  <  wscp  <  figm.  (A  good  approximation  to 
w lh  is 


which  is  always  less  than  Qgm.)  Qgm  will  also  play  an 
important  role  below. 

[15]  As  0  increases,  the  'P  curves  generally  shift  left  and 
down,  and  u)pcak  and  *P(ujpeak)  are  decreasing  functions  of  0 
as  confirmed  by  the  falling  curves  in  Figures  3a  and  3b. 
From  a  different  point  of  view,  if  for  some  value  of  0  the 
slope  3'P/3w  at  fixed  uj  is  positive,  the  slope  decreases  to 
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Figure  2.  The  function  'J'(w)  at  fixed  values  of  uipe/f^,  for  sin2  9  =  0.1  (top  curves),  0.3,  0.5,  0.7,  and 
0.9  (bottom  curves).  The  thick  solid  curves  in  Figures  2c  and  2d  are  for  0scp  and  have  zero  curvature  at 
u)  =  0.  Also  shown  in  Figure  2d  are  the  curves  for  Qgm,  which  crosses  zero  at  (upper  dashed 
curve),  and  for  0  =  ir/2,  which  crosses  zero  at  (lower  dashed  curve). 


negative  values  as  0  increases.  Extensive  numerical  tests 
show  that  once  this  slope  becomes  negative  it  stays  negative 
so  that  there  is  (at  most)  one  value  of  0  for  which  $  peaks  at 
any  fixed  u. 

[i6]  The  rising  curves  in  Figure  3a  show  the  inflection 
point  frequency,  where  $  switches  from  positive  curvature 
to  negative  curvature.  (The  values  wpeak  and  winn  are 
undefined  when  the  curvature  at  w  =  0  becomes  zero,  at 
sin2  0sep.)  This  ujjnfl  is  usually  an  increasing  function  of  0, 
although  it  must  stay  less  than  u)p<.ak.  In  fact,  Figure  3c 
indicates  that  Wjnfl  is  always  less  than 

4.1.  High-Density  Approximation 

[n]  In  this  and  the  next  subsection,  two  alternate  forms 
of  \F  are  presented  which  will  be  useful:  the  high-density 
approximation  and  a  novel  factorization  of  'F.  In  the  high- 
density  limit,  wpe/fi2  »  w/n„  simply  dropping  the  leading 
“  1  ’s”  in  the  definitions  of  R,  L,  and  P  leads  to  the  much 
simpler  expression 

^hd  _  M  f .  ^2  sin2  9 

uje  i  +  m  |  nen,-  2 

+  +  (i  -M)2cos2e.|  (is) 

Apart  from  the  terms  in  front  of  the  braces,  'iHD  is  the  form 
introduced  by  Lyons  [1974b]  and  used  for  quasi-linear 
diffusion  coefficient  calculations  almost  exclusively  ever 
since. 


[is]  It  is  proved  in  Appendix  A  that <  \Fwo.  Since 
admits  explicit  analytical  expressions  for  wpeak  and 
^(Wpcak),  't  provides  a  useful  upper  bound  on  the  true  \I> max- 

4.2.  A  Factorization  of  $ 

[19]  A  different  form  of  the  expression  for  $  is  suggested 
by  defining  w  =  u/fie,  and  R  =  R(b  and  similarly  for  L,  P,  S, 
and  D.  Then  the  full  f/,  with  the  same  choices  of  ±  as  for 
equation  (12),  can  be  written  as 


*-(t 


rL—PS  sin2  0  S 

- = - 1 f-  ~ 

D  2  P  D 


RL—PS\ 2  sin4  6 

D  J  4P2 


+  cos2  0 


(16) 


The  benefit  of  writing  \f<  in  this  form  is  that  for  u  >  f lgm,  the 
first  term  in  parentheses,  denoted  by  $+,  is  a  strictly 
increasing  function  of  w,  while  the  term  in  braces,  denoted 
by  \&_,  is  a  strictly  decreasing  function  of  u>.  A  proof  is 
outlined  by  Albert  [2004].  The  monotonicity  of  \P+  and  'F_ 
immediately  gives  the  bounds 


^+(wic)^-(wc/c)  <  'F(uj)  <  $+(u)yC)'P_(wiC),  (17) 


which  will  be  close  estimates  if  wLC  is  close  to  uuc.  The 
lower  bound  is  superfluous  because  the  minimum  value  of 
\I>  is  known  to  occur  at  either  or  w uc,  but  the  upper 
bound  is  valuable  when  wpcak  lies  between  w£C  and  uuc,  as 
discussed  in  section  5. 
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Figure  3.  (a)  The  frequencies,  normalized  by  Qe,  where  $  has  a  peak  (falling  curves)  or  inflection  point 
(rising  curves),  forwJe/ft2  =  10  (dotted  curves),  0.1  (solid  curves),  10-3  (dashed  curves),  versus  sin2  0. 
(b)  The  values  of  'I'  at  its  peak,  (c)  The  inflection  frequencies  normalized  by  wLH. 


[20]  It  is  worth  noting  that  the  heavy  ion  approximation 

M  =  0  reduces  ( RL  -  PS)/D  to  1,  which  simplifies 
equation  (16)  considerably.  It  is  proved  in  Appendix  A 
that  <  $  so  that 

*<"-°>  <  <  VHD.  (18) 

Also,  ¥f=0)  and  ¥^0)  are  monotonic  increasing  and 
decreasing,  respectively,  for  all  frequencies  0  <  w  <  min[fie, 

^pe]' 

[21]  Figure  4  shows  the  dependence  of  the  three  versions 
of  VP  on  w  for  12  different  combinations  of  upJVle  and  9.  In 
each  double  plot,  the  entire  range  0  <  w  <  Qe  is  shown,  as 
well  as  an  expanded  view  of  the  range  0  <  w  <  figm.  The 
heavy  dashed  line  is  the  full,  exact  TL  Also  shown  in  each 
plot  is  <SfHD  as  a  solid  line  that  always  lies  above  the  full  'P 
and  ¥M=0)  as  a  solid  line  that  always  lies  below  the  full  $ 
but  is  in  several  cases  almost  indistinguishable  from  it.  The 
vertical  scale  is  indicated  by  listing  the  peak  value  of  'VHD\ 
this  can  be  greater  than  1,  but  the  full  expression  for  $  is 
always  less  than  1. 

5.  Eliminating  Ranges  of  0 

[22]  If  the  integrals  in  equation  (5)  are  performed  over  w, 
the  resonance  condition  gives  an  easily  solved  quadratic 
equation  for  cos"  0  as  a  function  of  w.  Integrating  over  x  - 
tan  9  instead  involves  the  more  difficult  problem  of  solving 
for  w(9)  but  offers  greatly  improved  efficiency  by  cutting 
down  the  required  range  of  integration  or,  in  favorable 
cases,  eliminating  9  integrals  altogether.  As  in  Papers  1 
and  2,  the  strategy  is  to  circumvent  the  complicated  alge¬ 


braic  dependence  of  f  on  u  by  considering  'I'  geometri¬ 
cally,  using  the  properties  discussed  in  the  previous  section, 
while  exploiting  the  relatively  simple  algebraic  dependence 
of  $  on  9. 

5.1.  Overlap  of  Maximum  and  Minimum  Values 

[23]  The  fundamental  geometric  idea  is  that  at  fixed  0, 
F(w)  and  $(w)  cannot  intersect  if  Vis  larger  than  $  or  if  f7 is 
smaller  than  iP  for  all  w  between  wic  and  ujyC.  (The  latter 
case  is  illustrated  in  Figure  5a.)  Such  0  values  yield  no 
resonances  and  can  be  skipped  in  the  integration  for  the 
diffusion  coefficients.  The  prescription  is  thus 

Skip  0  :  Vmax  <  4* min  Or  Fmjn  >  'l;max  ■  (19) 

Once  the  values  of  w  that  maximize  and  minimize  V and  \P 
between  wiC  and  uuc  are  determined  for  fixed  0,  writing  out 
the  conditions  in  (19)  directly  gives  quadratic  inequalities  of 
the  form  A  cos4  0  +  5  cos2  9  +  C  >  0.  As  in  Paper  1 
(Appendix  Al)  or  Paper  2  (Appendix  B),  these  inequalities 
are  readily  solved  for  ranges  of  9  which  may  be  skipped. 
5.1*1.  Imijj  and  Fmax 

[24]  The  geometric  behavior  of  V,  shown  in  Figure  1, 
makes  it  easy  to  determine  Fmin  and  Fmax  (just  as  in  Papers 
1  and  2).  For  w„  =  0,  they  are  both  just  V0.  For  oj„  <  0,  Fmin 
is  F(wic)  and 

^max  is  F(u)jy c)-  For  0  <  <  Wyc>  ^min  IS 

V(uuc)  and  ^rnax  is  F(u>£c),  while  for  10 yc  <  Fmjn  is 
F(wiC)  and  Fmax  is  V(wuc).  If  wiC  <  H,  <  wuc,  Fmin  is  the 
smaller  of  V(^lc)  and  F(wyC)  and  Fmax  is  infinite. 

5.1.2.  The  Condition  Fmax<  'Pmin 

[25]  Since  ^mir  is  the  smaller  of  'P (ujic)  and  'P(oJyc),  the 
condition  Fmax  <  'Pmin  is  satisfied  by  9  ranges  where  both 
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Figure  4.  The  function  P(w)  (dashed  curves),  for  different  values  of  0  and  upJQ.e.  Also  shown  are  the 
high  density  approximation  (upper  solid  curves)  and  the  M  =  0  approximation  (lower  solid  curves).  For 
better  resolution,  the  ranges  0  to  Qgm  and  0  to  are  shown  separately.  The  vertical  scale  goes  from  0  to 
the  peak  value  of  PHZ>,  which  is  indicated  in  each  figure.  For  as  low  as  1,  the  high-density 

approximation  is  quite  good  for  small  w/fie  unless  0  is  very  large. 


Pmax  <  P(wiC)  and  Fmax  <  P(w[/C)-  Note  that  for  w„  <  0  or 
>  ujc/c,  Fmax  is  V(^uc)  and  decreases  with  increasing  |w|. 
Therefore  for  |«|  large  enough,  Fmax  will  always  be  less 
than  Pmjn,  which  gives  a  systematic  means  of  limiting  the 
doubly  infinite  number  of  n  values  in  the  sum  of 
equation  (5). 

5.1.3.  The  Condition  Fmin  >  Pmax 

[26]  This  condition  is  more  complicated,  because  Pmax 
depends  on  0.  From  Figure  3a,  uip(.ak  decreases  as  0  increases 
and  there  is  a  unique  value  0  =  Quc  (found  numerically)  for 
which  P  peaks  at  wye.  so  Pmax  =  P(u)t/c)  for  0  <  Quc- 
Similarly,  there  is  a  larger  value  0 LC  for  which  P  peaks  at 
uic,  so  Pmax  =  'P(uJic)  for  9  >  Qlc.  (Thus  wiC  <  ^uc  but 
©ic  >  Quc  )  For  either  9  <  Quc  or  0  >  Qlc  the  condition 
Fmin  >  'Pmax  specifies  quadratics  for  cos2  0. 

[27]  For  0  between  Quc  aud  Qlc,  the  maximum  of  P 
occurs  somewhere  between  wuc  and  u>lc,  and  no  simple, 
exact  characterization  has  been  found.  However,  any  func¬ 
tion  P  which  is  an  upper  bound  of  P  may  be  used,  since 
then  Fmin  >  Pmax  guarantees  Fmin  >  Pmax.  One  immediate 
(overestimate  is  P(wc/c,  Quc),  since  the  peak  value  of  P 
decreases  as  0  increases  (see  Figure  3b).  Note  that  for  9  = 
Qlc,  'Fpeak  is  actually  P+(wiC)  P -(wtc»  Qlc),  while  the 


estimate  just  given  is  P+(wj/c)  P_(wt/c>  Quc)  and  the  upper 
bound  given  in  equation  (17)  is  P+(wt/c)^-(wic.  Qlc)- 
Neither  estimate  dominates  the  other,  so  both  should  be 
used  to  narrow  down  the  retained  0  ranges.  The  upper  bound 
P™  can  be  used  as  well. 

5.2.  Curvature 

[28]  The  conditions  in  equation  (19),  based  on  the  max¬ 
imum  and  minimum  values  of  of  F(w)  and  P(w),  may  be 
augmented  by  geometrically  simple  ideas  involving  their 
curvature.  For  w  >  wlh,  P  has  negative  curvature  (concave 
downward)  and  F  has  positive  curvature  (concave  upward) 
on  either  side  of  (for  w„  positive)  or  for  w  <  |w„|/2  (for  io„ 
negative).  Now  suppose  V(mLc)  <  P(u>ic)  and  F(w uc)  < 
P(u)[yC)  for  some  range  of  9,  as  illustrated  in  Figure  5b. 
Since  F(w)  and  P(w)  “curve  away  from  each  other,”  there 
can  be  no  intersections  between  wLc  and  wye-  The 
corresponding  0  ranges  can  be  found  from  quadratics  as 
above. 

[29]  Similar  considerations  apply  to  the  combination  V 
concave  downward,  P  concave  upward,  with  F(wiC)  > 
P(u>ic)  and  V(u>uc)  >  PfoJ^c)-  However,  this  requires  both 

<  0,  uj  >  |w„|/2  (for  the  V curvature)  and  w  <  (for  the 
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Figure  5.  F(w)  and  'P(u)  with  to2e/fi2  =  1  and  sin2  9  =  0. 1 , 
for  a  0.5  MeV  electron  with  n  =  -1.  Also  shown  for 
illustration  are  uj£C  =  0. 1 8  and  -  0.36.  Maximum  and 
minimum  values  are  shown  by  filled  circles  for  V  and  open 
Circles  for (a)  a  =  85°.  In  this  case,  Fmax  <  tymin  and  there 
is  no  resonance  between  the  cutoffs,  (b)  a  =  80°.  Here,  Vmax 
is  not  less  than  but  the  lack  of  resonance  between  the 
cutoffs  may  be  detected  by  the  curvature  argument  of 
section  5.2. 


'P  curvature)  so  that  y  >  \n\Hy/M  or  E  >  10  MeV  for 
electrons. 

5.3.  Ranges  of  X 

[30]  The  diffusion  coefficients  may  be  bounce  averaged 
according  to  (D^ct0)  =  (1/t*)  /  Daa(d<V<9a)2<fr/v||,  etc., 
and  the  integration  along  a  field  line  may  be  carried  out  over 
latitude  X  [Lyons  et  al,  1972]  or  local  pitch  angle  a.  The 
preceding  approach  may  be  extended  to  eliminate  unneces¬ 
sary  .ranges  of  X  (or  a),  as  follows.  Consider  V  and  $  as 
functions  of  u>  =  u)/fie  so  that  the  curve  of  F(w)  at  fixed  9  is 
X-dependent  only  through  a  and  'P(w)  is  X-dependent  only 
through  u }pJ£le-  V  is  a  strictly  decreasing  function  of  X, 
while  the  X  dependence  of  \P  is  more  complicated  but  ’I'  = 
(ojpe/ff2)^  strictly  increases  with  upe/fie,  as  shown  in 
Figure  6  and  proved  in  Appendix  A.  Then  the  resonance 
condition  for  a  particle  with  equatorial  pitch  angle  Oo  can  be 
rewritten  slightly  as  V  =  ’F,  where 


V  & 


Jpe 0 


1  -  h  sin2  oo 


hi1 


(w  -  w„) 


:  COS2  9 


(20) 


and  h(\)  is  Lle  divided  by  its  equatorial  value  The 
bracketed  combination,  denoted  by  H(\),  strictly  decreases 
as  X  increases.  Any  X  dependence  of  oJpe,  w£C,  and  <^uc  will 
be  ignored  for  simplicity,  though  this  can  be  generalized. 
Then  the  normalized  frequency  cutoffs  and  wyc/ne 

are  bounded  below  and  above,  respectively,  by  w£c  = 
wzc/^e(^max)  and  &uc  =  (-Jj/c^e(^'min)>  where  Xmjn  and  Xmax 
are  determined  by  both  gq  and  the  range  of  the  waves. 

[31]  It  is  now  straightforward  Jo  derive  conditions  for  X 
ranges  which  guarantee  Fmax  <  I'min  or  Vmi„  >  $mia  so  that 
these  ranges  need  not  be  included  in  the  bounce  averaging. 
The  maximum  and  minimum  values  must  be  found  over 


ranges  of  9,_  X,  and  w.  For  simplicity,  consider  the  case  n  =  0, 
for  which  V  is  independent  of  u>.  An  underestimate  of  Vmm 
is  V  evaluated  at_9max  and_  Xmin,  and  an  overestimate  of 
I' max,  denoted  by  I' max.  is  ®  evaluated  at  9min  and  Xmin  and 
maximized  (numerically)  over  the  range  (bLc  to  il> uc-  Thus 
all  X  values  small  enough  that 

(v/c)2  (bJpeo/QeofHiX)  COS2  ^  >  $MAX  (21) 

will  satisfy  Vmin  >  Similarly,  all  X  values  large 

enough  that 

(v/c)2  (uvo/a<>)2//(X)  cos2  6mi„  <  Iw  (22) 

where  I'min  is  I'  evaluated  at  0max  and  Xmax  and  minimized 
over  uj£c  to  Cjuc,  will  satisfy  Kmax  <  ^min.  These  X  ranges 
will  have  no  resonances  within  the  imposed  9  and  w  ranges 
and  cannot  contribute  to  the  bounce-averaged  diffusion 
coefficients.  A  further  restriction  on  X  is  found  in  section  7. 

[32]  The  analysis  may  be  extended  to  n  ^  0  to  restrict  X 
ranges  and,  following  section  5.1.2,  to  identify  maximum 
values  of  |n|  beyond  which  there  are  no  resonances  at  any 
accessible  X. 

6.  Finding  Eligible  Resonances 

[33]  After  eliminating  as  much  of  the  9  and  X  integrals  as 
possible,  the  remaining  0  ranges  must  be  integrated  over 
according  to  equation  (5);  to  do  this,  the  resonances  that  do 
in  fact  lie  between  the  cutoff  frequencies  must  be  found.  A 
Straightforward  approach  is  to  derive  a  polynomial  P(w) 
from  the  resonance  condition  [Glauert  and  Home,  2005], 
which  turns  out  to  be  10th  order  in  w.  Instead,  the  geometric 
analysis  of  V  and  'I'  presented  above  may  be  used  to  try  to 
bracket  the  frequencies,  which  are  then  found  with  any  one¬ 
dimensional  root  finding  algorithm  [e.g.,  Press  et  al.,  1992]. 
The  simplest  situation  is  where  F'(w)  and  'F'(w)  have 


u/Qe  «/ ne 

Figure  6.  The  functions  (a)  ^  and  (b)  (w2e/ft2)  'P  versus 
u/Qe  for  different  values  of  w2e/fi2,  with  0  =  30°. 
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opposite  signs  everywhere  between  the  cutoffs,  since  then  F 
and  'I'  can  intersect  at  most  once  (and  will  do  so  only  if 
F(uiC)  -  ^(wic)  and  V{wug)  -  'P(wuc)  have  opposite 
signs).  Otherwise,  if  V  and  'I'  have  the  same  signs  but  V" 
and  \P"  have  opposite  signs,  then  V  and 'P'  can  have  at  most 
one  intersection,  which  separates  possible  intersections  of  F 
with  'I'.  These  cases  are  identified  below. 

6.1.  Case  it  -  0 

[34]  This  is  the  simplest  case  because  F  does  not  depend 
on  w.  For  0  >  0sep,  $(u)  is  always  decreasing  so  can 
intersect  F  only  once.  For  0  <  0sep,  there  will  be  at  most 
two  intersections  of  'F(lj)  with  V,  separated  by  wpeak.  In  fact, 
uipcak  will  only  lie  between  uiLC  and  uiuc  for  Quc  <  9  <  &lc- 
For  other  0  values,  ’P(ui)  is  always  either  increasing  or 
decreasing  between  the  cutoffs  and  intersects  F(uj)  there 
exactly  once  (the  0  values  for  which  these  intersections  are 
outside  the  cutoffs  have  already  been  eliminated). 

6.2.  Case  w„  >  0 

[35]  There  are  several  subcases.  First,  for  0  >  0?ep,  'I'  is 
always  decreasing,  so  it  can  intersect  the  rising  branch  of  F 
(where  w  <  oj„)  at  most  once.  An  intersection  with  the 
decreasing  branch  of  F  would  require  oj„  <  Qgm  so  that  y  > 
1  iVm  or  E  >  20  MeV  (for  electrons),  which  is  not 
considered  here. 

[36]  Otherwise,  $  is  either  rising  or  has  negative  curva¬ 
ture.  The  falling  branch  of  F  has  positive  curvature,  so  it 
could  intersect \P  at  most  twice,  with  the  two  intersections 
separated  by  the  single  root  of  V  -  IP'.  Of  course  these  roots 
need  only  be  sought  between  max(u)„,  u)ic)  and  uuc- 

[37]  Since  the  rising  branch  of  F also  has  positive  curva¬ 
ture  and  starts  out  below  >P  at  w  =  0,  it  can  intersect  'I'  at 
most  once  in  the  range  u  >  winfl  where  $  has  negative 
curvature.  However,  if  the  rising  branch  of  F intersects  \P  at 
oj  <  Wjnfl,  both  curves  are  rising  and  have  positive  curvature 
so  multiple  intersections  must  be  considered.  This  is  anal¬ 
ogous  to  the  situation  of  section  3.2. 1  of  Paper  1  and  is 
resolved  the  same  way:  V  is  necessarily  greater  than  'P'  at 
the  first  such  intersection  w*,  and  dV/du  increases  faster 
than  linearly,  so  even  the  function  (<9F/<9uj)/w  is  increasing 
with  u>.  However,  even  when  d'P/3u)  is  increasing,  (c?\P/&j)/w 
is  decreasing  (for  all  w,  verified  numerically).  Thus  V 
increases  faster  than  for  u  >  w.,  and  this  branch  of  V 
cannot  intersect  $  again.  In  summary,  there  are  at  most  three 
intersections:  one  for  the  rising  branch  of  Fand  two  for  the 
falling  branch. 

6.3.  Case  w„  <  0 

[38]  As  in  section  3.3  of  Paper  1,  this  is  the  most 
problematic  situation.  Cases  where  either  F(uj)  and  'P'(w), 
or  else  V  and  'P",  have  opposite  signs  everywhere  between 
the  cutoffs  can  be  handled  as  above.  However,  Fand  'P  both 
have  positive  slope  and  positive  curvature  for  w  <  min(w;nfl, 
|u>„|/2)  (so  that  uiiC  <  w £//),  and  they  both  have  positive 
slope  and  negative  curvature  for  max  (ijjinfl,  |oj„|/2)  <  u  < 
wpeak  (so  that  \n\  <  y  cos  0min).  If  the  preceding  tests  do  not 
rule  out  these  possibilities,  then  no  convenient  argument  is 
known  to  bracket  the  roots  of  F—  IP,  as  was  done  above  for 
w„  >  0  or  in  section  3.3  of  Paper  1. 

[39]  In  principle,  Sturm’s  theorem  gives  the  exact  number 
of  real  roots  of  the  polynomial  P( u>)  in  any  range  (ui|,  w2)  in 


terms  of  the  number  of  alternations  in  sign  of  sequences  of 
polynomials  derived  from  P(uj)  and  evaluated  at  U|  and  uj2 
[e.g.,  Barbeau,  2003],  However,  since  the  1 1  coefficients  of 
depend  on  the  parameters  wpe/Q.e,  0,  y,  a,  and  n  (as  do 
the  coefficients  of  the  derived  polynomials),  a  general 
characterization  of  these  signs  has  not  been  attempted. 
Descartes’  rule  of  signs  and  the  theorem  of  Fourier  and 
Budan  are  similar  but  less  definitive  and  almost  as  difficult 
to  apply. 

[40]  To  proceed,  one  option  is  to  revert  to  the  polynomial 
P(u)  for  o)(0),  accepting  the  inefficiency  that  the  solutions 
will  not  necessarily  lie  between  wLC  and  wuc-  Also,  the 
(real)  solutions  must  be  checked  to  verily  that  they  satisfy 
the  original  equation  F  =  'P.  An  alternative  is  to  integrate 
22“  over  u  instead,  finding  0(w)  easily  via  a  quadratic  in 
cos2  0  as  mentioned  in  section  5;  this  accepts  the  analogous 
inefficiency  that  much  of  the  w  range  will  likely  yield  0 
values  outside  the  range  0min  to  0max  and  sacrifices  the 
elimination  of  any  0  ranges  that  are  possible  for  this 
particular  value  of  uj„. 

7.  Dependence  of  Dpp  on  Density  Ratio 

[41]  The  ratio  of  momentum  diffusion  to  pitch  angle 
diffusion  is  controlled  by  equation  (4),  which  may  be 
written  as 


a  -  ^  =  cos2  a 


XT 

aa 


sin2  a  (w  -  f 1„f 


fi„  = 


(23) 


This  is  the  same  dependence  on  w  as  in  equation  (10)  for 
F(w),  with  the  similar  result  that  A  increases  with  increasing 
w  for  n„  <  0  or  w  <  Qn  and  decreases  with  increasing  u  for 
0  <  <  w. 

[42]  As  w pJSle  decreases,  Figure  6a  shows  that  \P 
increases  for  w  <  wpeak  while  the  F  curve  is  unchanged,  so 
(from  Figures  1  or  5a)  the  resonant  frequency  will  increase 
where  dV/du  >  0  and  vice  versa.  Thus  for  w„  <  0  or  w  <  ui„, 
decreasing  w pe/£le  means  typical  (<u>peak)  resonant  w  values 
increase  and  A  increases;  for  w  >  decreasing  w pJ^le 
means  typical  resonant  w  values  decrease  and  A  increases. 

[43]  The  only  remaining  case  is  0  <  <  w  <  f2„,  but  such 

resonances  do  not  seem  to  be  possible.  (Figure  5a  shows  two 
resonances  with  0  <  <  u>,  but  the  parameters  are  such  that 

w  >  n„.)  This  statement  is  proved  in  Appendix  B  with  the 
restriction  upe/Qe  >  1/2.  Extensive  numerical  tests  (looking 
for  valid  solutions  to  the  quadratic  for  sin2  0  for  over  101" 
combinations  of  w,  w pJVle,  7,  a,  and  n )  indicate  that  it 
holds  for  all  values  upeISle  >  10~3  (the  lowest  value  tested). 

[44]  This  finding  can  also  be  used  to  increase  the  effec¬ 
tiveness  of  the  methods  of  section  5.1,  since  for  w„  >  0,  only 
the  parts  of  the  interval  (wic,  wuc)  outside  (u>„,  fi„)  need  be 
considered.  This  finding  could  further  be  used  to  to  elim¬ 
inate  some  \  ranges,  as  in  section  5.3;  for  w„  >  0,  values  of 
X  for  which  the  interval  (ujiC,  uuc)  lies  within  (w„,  fi„)  can 
be  omitted. 


8.  Results 

[45]  Bounce-averaged  results  based  on  models  of 
CRRES  wave  data  covering  three  different  local  time 
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Figure  7.  Local  (equatorial)  diffusion  coefficients  Daa 
and  Dpp  for  a  1  MeV  electron  with  u )pJQe  =  1.5,  =  0.35, 

8w  =  0.15,  =  0.05,  and  ui[/c  =  0-65  (after  Home  et  al. 

[2003]).  The  solid  curves  are  based  on  the  full  whistler 
dispersion  relation  while  the  dashed  curves  use  the  high 
density  approximation. 

sectors  [Meredith  et  al.,  2001]  were  presented  by  Home  et 
al.  [2005],  The  diffusion  coefficients  were  calculated  by  two 
independent  codes,  one  described  by  Glauert  and  Home 
[2005]  and  one  based  on  the  analysis  presented  above.  In 
this  section,  the  reliability  of  the  high-density  approxima¬ 
tion  is  explored.  Also,  the  gain  in  computational  efficiency 
of  the  present  methods  is  demonstrated. 

8.1.  Comparison  to  the  High-Density  Approximation 
[46]  Local  diffusion  coefficients  were  calculated  for 
1  MeV  electrons  with  the  same  parameters  as  in  the  work 
of  Home  et  al.  [2003],  namely  S2(w),  a  gaussian  of  w/fle 
centered  at  wm  =  0.35  with  halfwidth  6w  =  0.15  and  cutoffs 
w LC  =  0.05  and  wuc  =  0.65,  and  gu(0),  a  gaussian  of  tan  0 


Figure  8.  Diffusion  coefficients  for  a  1  MeV  electron  as  in 
Figure  7  with  upJQe  =  1 ,  wm  =  0.4,  6w  =  0. 1 ,  wic  =  0.3,  and 
wt/c  =  0.5. 


Figure  9.  Diffusion  coefficients  for  a  1  MeV  electron  as  in 
Figure  7  with  wpJQ.e  =  0.5,  wm  =  0.1,  6w  =  0.05,  uiic=  0.05, 
and  W[/c  =0.15. 


centered  at  xm  =  0  with  halfwidth  xw  =  tan  30°  and  cutoffs 
*min  =  0  and  xmax  =  1.  For  upe!Vle  =  10,  7.5,  5,  and  2.5,  the 
values  obtained  for  Daa  and  Dpp  effectively  duplicate 
Figures  2  and  3  of  Home  et  al.  [2003].  Agreement  was 
also  obtained  for  upJVle  =  1.5,  although  this  required 
lowering  xmax  to  slightly  less  than  (i.e.,  0.999  times)  the 
resonance  cone  value.  Rather  surprisingly,  the  values 
obtained  using  the  form  of  the  diffusion  coefficients  based 
on  the  high  density  approximation,  computed  as  in  Paper  1, 
are  essentially  identical  for  wpe/tte  >  1.5  and  almost  the 
same  for  wpe/De  =  1 .5  (which  in  the  high  density  approx¬ 
imation  does  not  involve  the  resonance  cone),  as  shown  in 
Figure  7.  Figure  8  shows  another  comparison,  with  the 
normalized  frequency  parameters  changed  to  upe/Qe  =  1,  um 
=  0.4,  6w  =  0.1,  w ic  =  0.3,  and  < jjuc  =  0.5  (this  avoids  the 
resonance  cone);  the  high-density  approximation  holds  up 
quite  well  except  at  low  values  of  a.  Figure  9  shows  a 
comparison  with  =  0.5,  =  0.1,  6w  =  0.05,  u Jlc  = 

0.05,  and  wt/c=  015;  here,  the  high-density  approximation 
is  not  reliable. 

8.2.  Computational  Efficiency 
[47]  The  amount  of  computer  time  needed  to  perform 
these  calculations  varies  greatly  with  the  various  physical 
parameters,  as  well  as  numerical  details.  The  test  case  used 
here  is  the  calculation  of  bounce-averaged  diffusion  rates 
for  1  MeV  electrons  at  L  =  4.5  in  the  night  sector,  with  -5  < 
n  <  5,  as  presented  by  Horne  et  al.  [2005].  (Models  for  the 
other  sectors  considered  in  that  paper  approach  the  reso¬ 
nance  cone.)  The  CPU  times  reported  are  for  totals  for  ao  = 
89°,  88°,  . . .,  5°  (just  outside  the  loss  cone).  To  perform 
numerical  integrations,  the  code  described  here  relied  on 
Romberg’s  method,  adapted  from  the  implementation  of 
Press  et  al.  [1992],  The  calculations  were  performed  as 
/  /  dxD™,  where  x  =  tan  0.  As  a  baseline,  the 

resonant  frequencies  w(0)  were  found  by  solving  the 
polynomial  P(w)  discussed  in  section  6.3,  regardless  of 
n,  without  using  any  of  the  analysis  above.  As  expected, 
this  way  was  the  slowest.  Then  the  techniques  of  sections 
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Table  1. 

CPU  Times, 

s 

PM 

II 

E,  MeV 

All  0 

Restricted 

All  0 

Restricted 

0.01 

419.0 

M<5 

86.4 

147.4 

54.5 

0.03 

366.2 

91.0 

119.6 

51.3 

0.1 

364.7 

164.7 

127.3 

86.2 

0.3 

389.8 

263.4 

166.5 

140.5 

1.0 

508.3 

421.0 

235.1 

204.0 

0.01 

1248.5 

|n|  <  20 

86.4 

385.2 

54.5 

0.03 

1165.5 

90.4 

337.0 

51.3 

0.1 

1144.3 

165.5 

340.6 

87.3 

0.3 

1231.9 

331.7 

386.7 

172.4 

1.0 

1551.3 

834.5 

540.1 

396.3 

5.1  and  5.2  were  used  to  cut  down  the  ranges  of  0;  over 
the  ranges  not  skipped,  the  resonances  were  again  found 
using  P(uj).  This  produced  the  same  values  (within  a  few 
percent)  in  considerably  less  time.  Next,  solving  P(u)  for 
the  resonances  was  replaced  by  solving  V  -  9  but,  as  an 
experiment,  the  integral  over  wavenormal  angle  was  done 
over  the  whole  range  0  to  xmax.  This  also  sped  up  the 
calculations  (sometimes  by  more,  sometimes  by  less). 
Finally,  the  techniques  of  restricting  the  0  ranges  and 
solving  V-  9  were  combined,  yielding  the  fastest  method. 
The  timing  results  are  shown  in  Table  1.  Generally 
speaking,  lower-energy  particles  have  fewer  resonances, 
so  there  is  more  opportunity  to  avoid  nonresonant  values 
of  0,  n,  and  X.  For  the  same  reason,  the  techniques  will  be 
more  effective  the  smaller  the  range  of  allowed  resonan¬ 
ces,  to  ^uc- 

[48]  Because  of  the  decreased  computer  time,  it  is  feasible 
to  consider  larger  values  of  n.  For  large  E  and  small  oto,  it 
turns  out  that  resonances  are  present  with  n  as  large  as  ±23, 
although  the  addition  to  the  total  diffusion  coefficients  is 
only  a  few  percent.  Moreover,  the  method  discussed  in 
section  5. 1 .2  to  estimate  the  maximum  ]«|  needed  at  each  ao 
and  X  was  very  effective,  typically  overestimating  |nmaxl  by 
only  2  or  3.  Thus  the  code  takes  only  a  little  longer  to 
include,  in  effect,  — oo.  <  n  <  oo.  Because  the  simplest 
method  behaves  especially  badly  when  n  is  too  large 
for  resonance,  Table  1  shows  timing  comparisons  with 
-20  <  n  <  20. 

[49]  As  mentioned,  these  timings  are  highly  dependent  on 
the  parameters  of  the  problem  and  on  the  numerical  proce¬ 
dure.  However,  it  is  clear  that  the  analysis  and  techniques 
presented  here  can  greatly  increase  the  efficiency  and 
practicality  of  the  calculations  and  (for  a  given  amount  of 
computer  time)  improve  their  accuracy. 

9.  Summary 

[50]  For  whistler  mode  waves,  a  detailed  characterization 
was  given  of  $  =  1/p2  versus  w,  treating  0  as  a  parameter.  9 
was  shown  to  have  no  internal  minima  and  a  single  local 
maximum  whose  location  u>peak  decreases  monotonically 
with  0,  as  does  ^(ujpeak).  For  a  given  range  wLC  to  uuc,  this 
determines  the  value  of  'l/m,n,  as  well  as  the  value  of  0/maX 
when  >P(wpeak)  lies  outside  the  frequency  range.  Several 
estimates  were  given  for  9,^  when  ojpcait  lies  inside  the 


frequency  range,  including  one  based  on  the  factorization  of 
$  into  increasing  and  decreasing  functions  and  one  based 
on  the  high  density  approximation. 

[51]  The  minimum  and  maximum  values  of  the  function 
F(w)  were  also  found,  which  leads  to  conditions  on  0 
(equation  (19))  for  the  existence  of  solutions  to  the 
resonance  condition,  V  =  9.  This  and  related  constraints 
on  0  ranges  can  be  used  to  drastically  reduce  the  amount 
of  integration  over  0  required  to  evaluate  quasi-linear 
diffusion  coefficients.  A  maximum  value  of  |n|  can  be 
determined  for  which  there  are  any  resonances  at  all.  This 
approach  can  also  be  extended  to  identify  ranges  of 
latitude  X  which  contain  no  resonances;  this  was  carried 
out  in  detail  for  the  case  n  =  0. 

[52]  For  variable  and  parameter  ranges  that  do  contribute 
to  the  diffusion  rates,  with  >  0  it  was  possible  to  both 
enumerate  the  resonant  frequencies  and  to  isolate  them,  i.e., 
to  find  expressions  that  bracket  and  separate  them  so  that 
they  can  easily  be  found  numerically.  For  cases  with  <  0 
that  cannot  be  eliminated,  a  10th  order  polynomial  must  be 
solved  for  w(0)  or  else  the  integration  performed  over  u>, 
which  only  requires  solving  a  quadratic  for  sin2  0  as  a 
function  of  u.  These  techniques  were  used  to  evaluate 
diffusion  coefficients  for  relativistic  electrons  at  L  =  4.5 
using  realistic  models  of  the  wave  parameters,  as  reported 
by  Home  et  al.  [2005]. 

[53]  It  should  also  be  relatively  straightforward  to  modify 
this  analysis  for  other  types  of  plasma  waves.  In  particular, 
Z  modes  can  resonate  with  energetic  electrons  [Home  and 
Thome,  1998]  and  have  a  refractive  index  that  is  quite 
similar  to  that  of  the  whistler  mode  [e.g..  Carpenter  et  al., 
2003], 

[54]  A  brief  analysis  was  given  of  the  typically  increasing 
relative  strength  of  momentum  diffusion,  Dpp/D o^,  with 
decreasing  density  ratio  upeICle.  In  passing  it  was  observed 
that  for  w„  >  0  (n  <  0  for  electrons),  a  particle  with  pitch 
angle  a  has  no  resonances  between  u>„  and  u)„/sin2  a.  It  is 
tempting  to  speculate  that  this  is  somehow  related  to  the  gap 
in  chorus  wave  spectra  frequently  observed  at  Q.J2  [e.g., 
Anderson  and  Maeda,  1977]. 

[55]  Calculations  using  the  full  whistler  dispersion  rela¬ 
tion  were  compared  to  results  with  the  widely  used  high- 
density  approximation.  In  the  cases  chosen,  this 
approximation  was  found  to  be  quite  reliable  for  density 
ratio  w pJCle  >  1  and  qualitatively  good  at  iopJCle  =  1 
except  at  low  pitch  angle  but  failed  badly  at  upJVle  =  0.5. 
Finally,  the  computational  effectiveness  of  the  analysis  was 
demonstrated,  decreasing  the  computer  time  significantly 
for  calculations  with  —  5  <  n  <  5  and  even  more  so  for 
-20  <  n  <  20.  Since  the  analysis  is  able  to  accurately 
estimate  the  maximum  |n|  for  which  resonances  occur,  it  is 
possible  to  take,  in  effect,  —00  <  n  <  00. 

Appendix  A:  Proofs  of  Some  Properties  of  the 
Function  )P 

[56]  For  whistler  mode  waves,  with  P  <  0  and  R  >  0,  $ 
can  be  written  as 


Y sin2  0  +  S/RL  -  a\J Y2  sin4  0  +  (D/RL)2  cos 2  0,  (Al) 
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where  Y=  { RL  -  PS)/2PRL  and  a  is  the  sign  of  DL.  Taking 
the  (partial)  derivative  of with  respect  to  any  of  the 
parameters  of  {R,  L,  P)  gives 


=  Y'  sin2  9  4-  (S/RL)' 

YY'  sin4  0  +  <D/RL)(D/RL)'  cos2  0 
_  CT - -  -  ■  — 

y  Y2  sin4  0  +  (D/RL)2  cos2  0 


(A2) 


'I'  increases  as  p  decreases  and  w2/fl2  increases,  as  claimed 
in  section  5.3. 

Appendix  B:  Absence  of  Resonances  Between 
and  fl„ 

[6o]  With  upe/Qe  >  1/2,  it  is  readily  proved  that  there  are 
no  resonances  between  w„  =  |n|SVy  and  On  =  oj„/sin2  a. 
Since  $  <  9HD, 


This  may  be  rewritten  as 


pf  ft! 

4/'  —  —  sin2  0(1  -  cos <}>)  -  cos0sin<j) 


f  s  V 

sin2  0 

+  ytfj 

1 - - —  (1  -  cos  cj>)  -  cos  8  sin  <p 

(A3) 


where 


.  ,  ioY  sin2  0,  \D/RL\  cos  0} 

{cos<|>,  sin  <f>}  =  1  1 

Y  Y2  sin4  0  +  (.D/RL)2  cos2  0 


and  cos  0  sin  <])  >  0.  With  0  and  4>  considered  as  independent 
variables,  0  <  0,  <j>  <  it,  it  can  be  shown  that  the  bracketed 
term  in  (A3)  is  never  less  than  0.  The  form  (A3)  leads  to 
proofs  of  three  assertions  concerning  $  made  in  the  text,  as 
shown  below. 

Al.  Proof  That  >  ¥M=0) 

[57]  Taking  4>'  =  d'S/dM,  it  is  readily  checked  that  P1  and 
R'  are  negative  and  {S/RL)'  is  positive.  Therefore  d'Sf/dM  is 
never  negative,  which  proves  4/  >  4'(iU~0)  as  claimed  in 
section  4.2. 

A2.  Proof  That  *  <  VHD 

[58]  Next,  replace  the  leading  “1”  in  the  definitions  of 
R,  L,  and  P  with  an  additional  parameter  (3,  so  that  the 
exact  4>  is  given  by  (3  =  1  and  the  high  density 
approximation  i&HD  is  given  by  0  =  0.  With  4''  =  d^/d\ 3, 
(A3)  is  unchanged  but  now  P’  and  R'  are  positive.  Therefore 
4''  will  be  negative  if  ( S/RL )'  is  negative.  Evaluating 
(S/RL)'  gives  a  quadratic  function  f$>)  divided  by  a 
positive  (squared)  fiinction.  (A  computer  algebra  program, 
such  as  Mathematica  [Wolfram,  2003],  is  useful  here.) 
Both  the  leading  coefficient  and  the  discriminant  of  ftff) 
turn  out  to  be  negative  so  /is  downward-facing,  with  no 
real  roots.  Thus  /  is  always  negative,  {S/RL)'  is  negative, 
and  <94r/c/3  <  0,  which  proves  4/  <  ^HD  as  claimed  in 
section  4.1. 

A3.  Proof  That  (w^/fi2^  Increases  as  0 o2pJLl2e  Increases 

[59]  Finally,  consider  4/  =  (w2/n2)$.  Let  p  denote  fi2/w2e, 
and  define  (R,  L,  P)  as  (pi?,  pZ,,  p P).  Then  the  form  of 
4 >{R ,  L,  P)  is  exactly  the  same  as  that  of  4 '{R,  L,  P)  in 
(Al).  Furthermore,  since 

P  =  p-(1  + 

and  similarly  for  R  and  L,  here  p  plays  the  same  algebraic 
role  in  4'  as  0  did  in  4>.  The  calculation  of  <54'/<9p  <  0 
proceeds  exactly  like  that  of  94>/d0  <  0,  with  the  result  that 


4» 


max 


n2  cos2  0 

wpe  ^ 


<  cos2  0. 


(Bl) 


If  <  ne  then  l/-y2  <  sin4  a In2,  and  evaluating  Vat  u  =  Q„ 
gives 


Fmin  >  - - S‘n  COs2  ®  >  cos2  0-  (B2) 

1  -  sin  a 

If,  on  the  other  hand,  >  w„,  then  cos2  a  >  1  —  |n|/y 

and  evaluating  Kat  uj  =  Qe  (recall  'I'  <  0  for  larger  w)  gives 

Frnin  >  - - TV  COS2  0  >  COS2  0.  (B3) 

7  —  M 

In  either  case,  Vmm  >  4'max,  so  there  can  be  no  solutions 
to  V  =  4-. 
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